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Metapopulations are models of ecological systems, describing the interactions and the behavior of 
populations that live in fragmented habitats. In this paper, we present a model of metapopulations 
based on the multivolume simulation algorithm tau-DPP, a stochastic class of membrane systems, 
that we utilize to investigate the influence that different habitat topologies can have on the local and 
global dynamics of metapopulations. In particular, we focus our analysis on the migration rate of 
individuals among adjacent patches, and on their capability of colonizing the empty patches in the 
habitat. We compare the simulation results obtained for each habitat topology, and conclude the 
paper with some proposals for other research issues concerning metapopulations. 

1 Introduction 

The field of metapopulations ecology deals with the study of spatial systems describing the behavior of 
interacting populations that live in fragmented habitats [ 17 ]. The purpose of these models is to understand 
how the local and global dynamics of metapopulation systems, usually balanced between local extinc- 
tions and new colonizations of unoccupied patches, depend on the spatial arrangement of the habitat. 
Consequently, relevant insights into related fields of ecological research, such as evolutionary ecology or 
conservation and landscape management, can be achieved. Indeed, the topology of fragmented habitats 
potentially holds relevant implications for the persistence of populations, and their robustness against 
natural or anthropogenic disturbance [36]. 

Recently, in addition to ever increasing applications of graph-based methods for the analysis of com- 
plex networks in cell biology [T] [2j , graph theory has also been applied to the study of metapopulations 
systems. In graph models of metapopulations, nodes are used to represent habitat patches, and graph 
edges are used to denote some functional connections between patches (typically related to the dispersal 
of individuals). Attributes can be associated to nodes, describing the quality or dimension of patches, 
while different types of edges can be exploited to represent the distance between connected patches, the 
rate of dispersal between a couple of patches, or simply whether two patches are connected or not. 

Metapopulation models using graph-based methods |[36l [T31 are simple to implement and require 
relatively few data for their definition, while individual-based models implement more detailed aspects, 
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concerning the nature and the interaction of populations ||34l l4l. Both types of modeling approaches are 
useful for the analysis of specific features of metapopulations but, while the first focuses on the properties 
of the habitat topology, the second is more concerned with the emergent dynamics. In this paper, we 
present a stochastic multi volume model of metapopulations, which integrates the explicit representation 
of interactions between the individuals of the populations - and therefore allows to simulate the emergent 
local and global dynamics - with a graph description of the habitat topology - which allows to investigate 
the influence of distinct spatial structures on the dynamics. 

This model, which represents a simplified extension of a previous metapopulation model that we 
introduced in Q |6l, is based on the multivolume stochastic simulation algorithm tau-DPP ifTTl [H, a 
stochastic class of membrane systems. Membrane systems, or P systems, were introduced in ll27l as a 
class of unconventional computing devices of distributed, parallel and nondeterministic type, inspired 
by the compartmental structure and the functioning of living cells. The basic model consists of a mem- 
brane structure where multisets of objects evolve according to given evolution rules. A comprehensive 
overview of P systems and of its many applications in various research areas, ranging from Biology to 
Linguistics to Computer Science, can be found in |[28l[T2ll29l . 

In tau-DPP, the distinct compartments of any multivolume model can be arranged according to a 
specified hierarchy (e.g., a membrane structure), under the additional assumption that the topological 
structure and the volume dimensions do not change during the system evolution (each volume is assumed 
to satisfy the standard requirements of the classical stochastic simulation algorithm, see lfl6l and Q for 
more details). Inside each volume, two different types of rules can be defined: the internal rules, which 
modify the objects contained inside the volume where they take place (in the case of metapopulation, 
they describe the growth and death of population individuals according to the Lotka-Volterra model of 
preys and predators), and the communication rules, which are used to move the objects between adjacent 
volumes (in the case of metapopulation, they describe the migration of population individuals). 

In this paper, tau-DPP is exploited to analyze the emergent dynamics of metapopulation systems, 
where the focus is on the influence that the topology of patches has on the migration of individuals, and 
their capability to colonize other patches in the habitat. To this purpose, we consider six different habitat 
topologies, formally described by graph structures, and analyze how the topological structure of patch- 
to-patch connections, and the rate of individual dispersal between connected patches, influence the local 
and global dynamics of a metapopulation. In particular, we will first consider how a given topology and 
a fixed dispersal rate between patches can influence the prey-predators dynamics, and then we will focus 
on the colonization of empty patches, starting from the dispersal of predators that live in a few patches 
which occupy peculiar positions in the given network topology. 

The paper is structured as follows: in Section [2]we present the concept of metapopulations in Ecol- 
ogy, and then describe the multivolume model of metapopulations by focusing, in particular, to the 
different habitat topologies. In Section[3]we will show the simulation results concerning the influence of 
these habitat topologies on the emergent dynamics of metapopulations, considering the effects of preda- 
tors dispersal and colonization. Finally, in Section [4] we conclude the paper with some final remarks and 
several proposals for further research issues concerning metapopulations. 

2 Metapopulations 

In this section, we first provide a brief introduction to the most relevant features of metapopulations, 
concerning both the topology of the habitats and the emergent dynamics. Then, we describe the modeling 
approach used in this paper, that is based on a stochastic class of membrane systems, which will be used 
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in Section|3]to analyze the influence of different network topologies on the dynamics of metapopulations. 
2.1 Dynamical models of interacting populations in Ecology 

Since its introduction in Il22ll . the concept of metapopulations (also called multi-patch systems) has been 
extensively applied in Ecology to analyze the behavior of interacting populations, to the purpose of deter- 
mining how fragmented habitats can influence various aspects of these systems, such as local and global 
population persistence, or the evolution of species lfl"8ll . Lately, this topic has been largely employed for 
other populations species, living in both natural and artificial/theoretical fragmented landscapes lPT7ll . 

A metapopulation consists of local populations, living in spatially separated habitats called patches 
- which can be characterized by different areas, quality or isolation - connected each other through a 
dispersal pool, which is the spatial place where individuals from a population spend some lifetime during 
the migration among patches. In multi-patch systems, two principal types of dynamics exist: on the one 
hand, the individuals of the different populations can have local interactions inside each patch (according 
to a given dynamical model, e.g., the Lotka-Volterra system of interaction between preys and predators 
ll25l ): on the other hand, the dispersal of individuals among mutually connected patches can influence 
the global behavior of the whole system II201I2T1I331I371 . The dispersal of individuals, which is usually 
dependent on the distance between patches, may reduce the local population growth, and thus increase 
the extinction risk, which can be due also to environmental and demographical stochasticity. Hence, 
the persistence of populations is assumed to be balanced between local extinctions and the process of 
colonization, that is, the establishment of new populations in empty patches |[T7l . 

Several theoretical frameworks for metapopulation analysis have been defined up to now, remarking 
specific properties of multi-patch systems which have been either explicitly or implicitly considered in 
these modeling methods (see, e.g., lfl4l [171 l24l [191 for further details). For instance, referring to the 
landscape, most theoretical models take care of the spatial structure of the habitat, the local quality of 
the environment, the patch areas and their mutual connectivity (or isolation), in order to capture the 
effect of habitat fragmentation on species persistence. In fact, good local conditions can determine the 
growth and the survival of populations inside the patches, and high patch connectivity can decrease local 
extinction risk. Moreover, as dispersal and colonization are distance-dependent elements, they can be 
used to account for the importance of real landscape structures. Referring to population interactions and 
dynamics, colonization can depend or not on the cooperation of migrating individuals (in the first case, it 
is called "Allee effect"). Models not accounting for within-patch dynamics - but only assuming whether 
a patch is occupied or not - usually consider local dynamics on a faster time scale with respect to the 
global dynamics, and also neglect the dependence of colonization and extinction rates on population 
sizes. Finally, regional stochasticity can account for "bad" or "good" years over the local environmental 
quality, which depends on, e.g., the weather conditions which affect sustenance resource availability and, 
once more, they can influence the growth and survival of populations. 

Recently, graph-based models for metapopulations have started to be more and more defined be- 
cause of the intuitive and visual way they hold for the representation of these ecological systems (see 
ll36l l23l l35l and references therein). In these models, nodes represent habitat patches and graph edges 
denote functional connections between patches (typically related to the dispersal of individuals). In ad- 
dition, attributes can be associated to nodes, describing the quality or dimension of patches, and different 
types of edges can be adopted to represent the distance between connected patches, the rate of dispersal 
between a couple of patches, or simply whether two patches are connected or not. These models allow 
to make insights into the features of habitat distribution, such as the predominant importance of some 
nodes or clusters of nodes with respect to other characteristics of metapopulation, like their dynamics, the 
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vulnerability to disturbance, the persistence of populations according to dispersal, and so on. These re- 
sults open promising perspective in related research fields as evolutionary ecology, conservation biology, 
epidemiology, management and design of natural reserves. 

2.2 A P system-based model of metapopulations: focusing on network topologies 

Most of the issues discussed in Section lZTl were explicitly considered in our previous model for metapop- 
ulations 13. In those works, metapopulation models were based on a class of membrane systems 
called DPP |[3T1 l30l . which were used to execute qualitative stochastic simulations of the local and 
global dynamics of metapopulations. In particular, in Q we introduced a model of metapopulations 
with predator-prey dynamics, where additional features were used in order to catch and better describe 
relevant properties of the modeled system. For instance, the regions of the membrane structure were 
represented as nodes of a weighted graph with attributes, where the weight associated to edges cor- 
responds to the "distance" among connected regions, while attributes specify their surface dimension. 
These new features are necessary in order to outline the spatial distribution of patches and the relevant 
additional features associated to them: the dimension of a patch is needed to define the density of the 
populations living inside that patch, while the distance is needed to identify isolated patches, as well as to 
define the dispersal rates of migrating individuals. Moreover, by using some rules which do not modify 
the objects on which they act (the so-called "mute rules"), we modified the classical view of maximal 
parallelism, by allowing the maximal application of rules but, at the same time, reducing the maximal 
consumption of objects. The model was applied to investigate some emergent metapopulation behaviors, 
such as the influence of patch dimension, patch-to-patch distance, stochastic breeding, the dynamics un- 
derlying migration and colonization, the effects due to isolated patches, etc. Then, in 01 we extended 
the analysis of that model by focusing on periodic resource feeding strategies, and compared different 
systems where either increasing, decreasing, stationary or purely feeding stochastic phases were defined 
inside each patch. We have shown there, for instance, how the seasonal variance can transform the basic 
Lotka-Volterra dynamics inside each patch into a more complex dynamics, where the different phases of 
a feeding cycle can be identified through the effect that they have on the standard oscillations of preys 
and predators. 

In this section, we present a simplified model of metapopulations, which exploits the multivolume 
stochastic simulation algorithm tau-DPP |[TTl l5l. With respect to the previous model, here we will not 
need to use the concept of mute rules, as the probabilistic choice and applications of rules is already 
embedded in the tau leaping algorithm |[T0l . on which tau-DPP is based. Moreover, we will not consider 
the presence of the dispersal pool, but we will instead focus our analysis on the direct communication 
of individuals among interconnected patches, according to some fixed network topologies. In order to 
compare the influence of each network, we have decided to perform our analysis on a total of 6 patches, 
spatially arranged in different ways. Namely, we assume that these network topologies can be described 
by graphs having the same number of nodes, but distinct connections, such as the chain, grid, star, ring, 
complete or random structure (see graphs a,b,c,d,e,f, respectively, in Fig. [[])■ From now on, we will 
refer to the formal data structure by using the term 'graph', and use the term 'network' to denote the 
topological relationship on each graph. 

Formally, each network topology v G {a,b,c,d,e,f}, can be generally described by a weighted 
undirected graph G v = (N^,E V ,w v ) where: 

• N£ is the set of nodes, such that each node pi G N£, i=l, . . ., 6, is characterized by a value 8 (pi) G A 
(with A being a set of attributes of some kind); 

• E v C {(pi,Pj) | PhPj G is the set of (undirected) edges between nodes; 
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(a) chain (b) grid (c) star 




(d) ring (e) complete (f) random 




Figure 1: Network topologies. 
• w v : E v — >■ IR + is the weight function associating a cost to each edge. 

In the case of metapopulations, the set of nodes N% coincides with the set of patches, the attribute 
of a node represents the area of the patch, the edges characterize which patches are directly reachable 
from any patch (self-edges might exist as well but will not be considered in this work), and the weight 
wjj of an edge (pi,Pj) represents a cost to measure the effort that individuals have to face when moving 
from patch pi to pj. Given a network topology v, we denote by Adj{pt) v the set of nodes that are 
directly connected to any node pi, that is, Adj(pi) v = {pj G N£ \ 3 (pi,Pj) G E v }. We also denote 
by deg(pi) v the degree of patch p,-, that is, the number of patches directly connected to pi (formally, 
deg(pi) v = card{Adj{pi) v )). We outline that, in what follows, we will assume that: (1) wjj = 1 for 
each (pi,pj) G E v and each v G {a,b,c,d,e,f}, that is, all edges have the same cost; (2) 5(p,-) = 1 for 
each pi G N£ and each v G {a,b,c,d,e,f}, that is, all patches have the same dimension. The rational 
behind this is that, in this paper, we focus our attention on the influence that different topologies of the 
habitat network can have on the local and global dynamics of metapopulations, regardless of the local 
features of each patch, or of the distances between patches. These features might be naturally added 
in further works related to this model, where real data can be used to define a specific model of some 
metapopulation systems. 

In addition to the chosen network topology, this model of metapopulations also considers the pres- 
ence of species individuals, which locally interact according to a chosen dynamics, and give rise to global 
dynamics thanks to the dispersal processes. To this purpose, in this paper we assume that each patch is 
characterized by the Lotka-Volterra (LV) model describing the interaction between the individuals of two 
populations, namely preys and predators. Inside each patch, the LV model is described by the following 
set of internal rules: 

n : AX -> XX 
r 2 : XY -> YY 
r 3 : Y -> X 

where X denotes the preys, Y denotes the predators, A denotes the sustenance resources and A is the 
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empty symbol. Rules r\ and r 2 model the growth of preys and predators, respectively, while rule r 3 
models the death of predators. Each rule is also characterized by a stochastic constants (expressed in 
time -1 ), that is used - together with the current amounts of individuals occurring in the patch - to evaluate 
its application probability step by step, according to the tau leaping algorithm (see ifTUl [TT1 HI for more 
details). All the simulations shown hereafter have been executed using the following values of stochastic 
constants and of initial amount of preys, predators, and sustenance resources: c\=0. 1, C2=0.01, C3=10, 
Xo=Fo=1000, Ao=200 (the value of A is fixed for the entire duration of each simulation). The simulations 
have been performed with the software BioSimWare (H, that implements different stochastic simulation 
algorithms for both single and multivolume systems. The software is available for free download at 
http://bimib.disco.unimib.it/index.php/Software. 

In Fig. |2] we show the oscillating dynamics (left side) of preys and predators in the single patch, 
obtained with this choice of parameters, and the corresponding phase space (right side). These figures 
can be considered as reference to compare and discuss the dynamics obtained in the multi-patch model, 
as described in Section [3] 
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Figure 2: The Lotka-Volterra dynamics in the single patch: oscillations in preys, X, and predators, Y (left 
side), and corresponding phase space (right side). 

The single patch model is then extended to a multi-patch model where, inside each patch /?, of each 
network topology v, we add as many communication rules as the number of patches connected to pi 
(that is, a total of deg{pi) v rules inside each patch). These rules are needed to move population indi- 
viduals among the various patches of the network, thus allowing to analyze the effects of migration and 
colonization in the metapopulation. This is done by attaching a destination target to each communication 
rule, specifying the destination patch, as it is usually done in P systems. Formally, in each patch pi of 
network v, we add the so-called dispersal rules 

rd Pj : Y — >• (Y, target (pj)), 

for each pj € Adj(pi) v . Similarly to the local rules ri,r2,rj, the probability of applying each dispersal 
rule is determined by using its stochastic constant Cd p , whose values will be given in the next section to 
consider different migration rates. 
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3 The influence of network topologies on metapopulation dynamics 

In this section we analyze how the topological structure of patch-to-patch connections, and the rate of 
individual dispersal between connected patches, influence the local and global dynamics of a metapop- 
ulation. In particular, in Section 13.11 we consider how a given topology and a fixed dispersal rate can 
influence the prey-predators dynamics, while in Section [3721 we focus on the capability of colonization 
of empty patches, starting from the dispersal of predators living in a few patches which occupy peculiar 
positions in the given network topology. 

3.1 Network topologies and migration 

In this section, we analyze the role of migration and compare the six network topologies with respect to 
four different conditions for the dispersal rules. Namely, we assume that each patch of each topology is 
initialized with a complete LV model as given in Section |Z21 where the value of the stochastic constant 
Cd p . for the dispersal of predators, in each patch pi £ N%, can assume one of the following values: 

1 . c dp . = 1 , for each pj G Adj (pi) v ; 

2. c dpj = 1 0, for each pj e Adj(p t ) v ; 

3. c dp .=20, for each pj eAdj(pi) v ; 

4 " %=3TO ' for each PJ G Ad J(Pi) V - 

By considering the first condition as reference, the power of dispersal in the second (third) condition is 
ten-fold (twenty-fold) the first one, irrespective of the position that patch pi occupies in the considered 
network. In other terms, the flux of dispersal from each patch, in the first three conditions, results 
amplified by the number of connections that each patch has with respect to the other patches in the 
network. On the contrary, the fourth condition corresponds to the situation when, for each patch pj € 
Adj(pi) v , the sum of the values of constants of dispersal rules in pi is always equal to 10, but the rate 
of dispersal along each edge from p\ to pj depends on the degree of p\. For instance, in the network 
topology a (Fig. [B, the value of ca. in patches po and p$ is equal to 10, while the value of c d in 
patches p\, . . ., p4 is equal to 5; in the network topology c (Fig. [T]), the value of cj in patch po is equal 
to 2, while the value of in all other patches is equal to 10, and so on. So doing, we can weigh the 
dispersal of predators according to the position of each patch in the network, and simulate a situation 
where the flux of dispersal from each patch towards its adjacent patches is uniform throughout the whole 
network. 

For space limits, in Fig. [3] we present the phase spaces of all network topologies, obtained from 
simulations of the fourth condition only. For each network, in particular, we show the phase space of the 
local dynamics of each patch. The graphics show that, in the case of the chain graph (phase space (a)), the 
patches having different degrees are characterized by different dynamics: in fact, patches po and p$ show 
a different behavior with respect to the other patches. In addition to the role of patch degree, we can see 
that also the position of patches in the graph plays a central role: despite the fact that patches p\,p2,P3 
and P4 have all the same degree, the dynamics inside p\ and p 4 differs from that of patches p 2 and p^. 
This is due to the different power of dispersal rules of their two neighbors, namely c</ = 10 in patches 
po, ps, while en = 5 in patches p 2 , ps, which cause a larger flux of predators dispersal towards patches 
p\ and p4. The global effect is the presence of three different dynamics (one in po, p$, another one in p\, 
P4, and a third one in p2, p-$), all of which are characterized by oscillations in X and Y with no regular 
amplitudes (compare these phase spaces with the standard LV phase space in the single patch model 
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Figure 3: The power of migration: LV dynamics in the phase space of each network topology. 
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given in Fig. |2 right side, and also with the phase spaces in Fig. [3j graphics (d) and (e)). Furthermore, 
we can evidence that these oscillations are characterized by an initial wider amplitude, which is reduced 
during time. 

Similarly, the dynamics of the patches in the grid graph (phase space (b)) is influenced only by the 
number of edges; in this phase space, we can identify two different types of dynamics: one for the patches 
with three edges (pi, p$) and another one for those with two connections. 

In the star graph (phase space (c)), the LV dynamics endures in all patches apart from po, where the 
number of preys X collapses to an attractor in zero, and no oscillations according to the LV dynamics 
in both X and Y can be established. In this patch, the number of predators fluctuates in a certain range, 
because of their dispersal from/to the other patches. Basically, in this condition patch po, that represents 
the center of the star, becomes a local area of the habitat where only dispersal occurs. 

The simulations for the ring and complete graphs (phase spaces (d), (e)) show very similar results: 
in both cases, all patches in each graph have the same degree (two in the first configuration and five in 
the second one), leading to regular oscillations in X and Y with almost constant amplitude. 

The results concerning the last configuration, the random graph (phase space (f)), show a combination 
of the effects described above. In particular, the dynamics of the patches differ each other depending on 
the degree of the patches themselves; moreover, in p4, which is characterized by the highest degree, the 
high number of incoming predators (migrating from the four adjacent patches) leads to the extinction of 
preys (similarly to what happens in patch po of the star graph). 

We also tested, for each network topology, the other three conditions listed above. In these cases, 
the results have shown that the amplification of the power of dispersal with respect to the patch degree 
gives rise to a balance between the incoming and migrating individuals, which leads to comparable LV 
dynamics for all networks, with regular oscillations inside each patch (data not shown). 

3.2 Network topologies and colonization 

In this section, we compare the six network topologies with respect to the capability of colonizing the 
empty patches that each network contains, starting from the patches that contain a complete LV model 
and that occupy a peculiar position in that network . We recall that in this work we are considering only 
the migration of predators, hence the empty patches are hereby assumed to contain no predators but only 
an initial amount of preys. In each network v, the set of patches initialized with the complete LV model 
will be denoted as p v LV . To test the feature of colonization, we consider four different initial conditions, 
hereby denoted as ICk, k=l, . . . ,4, where Yq=0 and: 

1. IC1 is characterized by c</ p =1 andXo=10; 

2. IC2 is characterized by c^.=l and Xo=100; 

3. IC3 is characterized by =10 and Xo=\0; 

4. IC4 is characterized by Cd p =10 and Xo=100. 

In each given network, all empty patches are initialized with the same chosen condition ICk, besides the 
patches in the set p v LV that are initialized with a standard LV model, having the communication constant 
Cd p . equal to the one given in the chosen ICk, and all other parameters as given in Section I2T21 

With this type of analysis, we expect to determine which features of the network topologies are more 
relevant with respect to the colonization of empty patches, under a given initial condition. All conditions 
have been tested for each network and, for each fixed initial condition, different sets of p v LV have been 
considered. In the following, for space limits, we present only some results of these simulations, and 
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briefly discuss the results obtained in the other analyzed conditions. In each of the following graph, 
preys (X) are represented with solid lines, while predators (Y) are represented with dashed lines. 

We start by considering the network v = a, that is, the chain graph. In this case, we present the 
results obtained in all the initial conditions IC1, IC2, IC3, IC4, considering three sets of LV patches, 

namely Plv={po,Ps}, Plv={pi\ and Plv={po}- In the fi rst case (plv={Po,Ps}, shown in Fig. H) we 
can see that, when the power of dispersal is low (IC1, IC2), the time required by the predators to reach 
the patches p2 and p^, which are at the highest distance from po and p^, allows an initial uncontrolled 
growth of the preys in P2 and p^, which subsequently undergo extinction as soon as the predators enter 
the patch. Such "delay" in the local establishment of a population of predators is the effect that prevent 
the formation of the LV dynamics; this effect, as shown hereafter, is a common aspect of all network 
topologies. Concerning the chain network, this is more evident in condition IC2, where the initial amount 
of preys inside the empty patches is higher than IC1: in this case, the LV dynamics can be established 
only in four of the six patches. On the other hand, with the initial conditions IC3 and IC4, the power 
of dispersal is sufficient to colonize all of the patches, irrespectively of the numbers of preys that are 
initially present in the empty patches and of the position of the LV complete patch. Similar results for 
the chain network have been obtained in the second analyzed case {Piv={pi\, shown in Fig. [5]> and in 
the third case {p c [y={po), data not shown). 





Figure 4: Colonization in the chain topology, with p" lv ={pq,Ps} and initial conditions IC1 (top left), IC2 
(top right), IC3 (bottom left), IC4 (bottom right). 
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Figure 5: Colonization in the chain topology, with pfy={p 2 } and initial conditions IC1 (top left), IC2 
(top right), IC3 (bottom left), IC4 (bottom right). 



For the network topology v = b, that is, the grid graph, we show the results obtained in the cases IC1, 
when Pi V ={po} (Fig. [6l left side) and p l [ v ={p l } (Fig. [6l right side). According to the position of the LV 
complete patches in this network topology, we can see that, in the first case, the predators are capable to 
colonize patches p\ and p?,, that are directly connected to po, and patch P4, that is directly connected to 
both p\ and p^. However, patches p 2 and ps cannot be colonized. In the second case, the higher degree 
of the LV complete patch p\, allows the colonization of all patches. With the initial condition IC2 (data 
not shown), in the other tested cases p l [ v ={p } and p h LV ={pi}, only the patches directly connected to po 
and pi, respectively, are colonized by the predators. 

For the network topology v = c, that is, the star graph, we show the results obtained in the cases 
IC1, when pjy={pi} (Fig. |7J left side) and p c LV ={p\,p^} (Fig. |7J right side). According to the position 
of the LV complete patches in this network topology, we can see that, in the first case, no patches are 
colonized because of the high degree of po (which is the only patch connected to p\) that spreads the 
predators over the other patches, thus preventing the formation of the LV dynamics. In the second case, 
the combined effect of migration from p\ and jt>3 allows the colonization of patch po, which is directly 
connected with both of them. We then performed other simulations starting with conditions IC3 and 
IC4: in these cases, the higher value of Cd v allows the colonization of every patch (except from patch po) 
independently from the initial position of the LV complete patch (data not shown). On the contrary, when 
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Figure 6: Colonization in the grid topology, with initial condition IC1 and Pi V ={po} (left), P b Lv = {Pi\ 
(right). 





Figure 7: Colonization in the star topology, with initial condition IC1 and p c LV ={pi} (left), pjy={pi,P3} 
(right). 





Figure 8: Colonization in the ring topology, with pfy = {Po} an d initial condition IC1 (left) and IC2 
(right). 
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we assume p c LV ={po}, that is, the center of the star, then all patches are fully colonized, independently 
from the considered initial condition. 

For the network topology v = d, that is, the circular graph, we show the results obtained in the cases 
IC1 and IC2, when p c [ v ={p } (Fig. [U left and right sides, respectively). Starting with the initial condition 
IC2, the predators are capable of colonizing only the patches directly connected to the LV complete patch 
Po, while in the case IC1, also patch P4 (being at distance 2 from the LV complete patch) is colonized. 
These results highlight, in particular, another aspect that was more marginal in the other simulations: the 
stochastic nature of the communication process and of the growth of preys, which leads to the extinction 
of preys in patch p2, while in patch P4 it drives the local behavior to an oscillatory dynamics. 

For the network topology v = e, that is, the complete graph, we show the results obtained in the cases 
IC1, when p e LV ={p } (Fig. |9l left side) and p e LV ={po,P3} (Fig. |9l right side). While in the second case 
- where the LV dynamics is initially placed in two patches - the predators can colonize all patches, in 
the first case the colonization of all empty patches fails. Once more, this is an effect of the stochastic 
noise combined with the low amounts of predators, which is in turn caused by the fact that the higher the 
number of adjacent patches, the lower the number of predators that persist inside each patch. In all other 
simulations performed with initial conditions IC3 and IC4, all patches have always been colonized, as 
the higher values of dispersal rules assure a more uniform spread of predators throughout the network, 
and thus flattens the influence of migration delay (data not shown). 

For the network topology V = f, that is, the random graph, we show the results obtained in the cases 
IC1, when pjy={p } (Fig. [lOl left side) and p^ v ={p 2 } (Fig. \\0\ right side). According to the position 
of the LV complete patches in this network topology, we can see that, in the first case, all patches are 
colonized by predators (similar results are obtained by placing the LV complete model in patch P4 - data 
not shown). In the second case, patch p$ is not colonized because there is only one path of length 2 which 
connects it to the initial complete LV patch p2\ the same holds for patch pi, which has distance from p 2 
equal to 3. For similar reasons, considering the case of initial condition IC1, with the LV complete model 
in patch p?,, the only patch that is not colonized by predators is p 2 (data not shown). In all the simulations 
performed with the initial condition IC2, some of the patches have not been colonized because of the high 
amount of preys initially occurring in the patches. On the other hand, with the initial conditions IC3, IC4, 
the power of dispersal allows the colonization of all patches (data not shown). 




Figure 9: Colonization in the complete topology, with initial condition IC1 and p e LV ={po} (left), 
P e LV ={Po,P3} (right). 
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Figure 10: Colonization in the random topology, with initial condition IC1 and p LV ={po} (left), 
Plv={P2} (right). 

4 Discussion 

The fragmented habitats of real metapopulations are usually characterized by complex network topolo- 
gies. In this paper, we have analyzed six small topologies that can be considered representative of local 
areas in a structured habitat, and we have investigated the influence that the degree and the position of 
each patch in the topology can have on the migration of individuals, as well as on the capability of col- 
onizing empty patches. Our analysis suggests that, with respect to the power of migration (Section [37Tb . 
we can identify different behaviours that depend on two characteristics of the topology: on a first level, 
the local behaviour inside each patch is influenced by its degree. This is especially evident if we compare 
the network topology described by the circular or complete graphs, with the topology described by the 
star graph: while in the first case (where all nodes have the same degree) all patches are characterized by 
a similar (regular) oscillating dynamics, in the second case the most critical node is the center of the star 
(which has a much higher degree than all other nodes in the same graph). In the latter case, this patch is 
likely to undergo a local modification of its initial dynamics, due to a more higher incoming migration 
of individuals from all other adjacent patches. On a second level, assuming in this case that the degree 
of nodes is equal, then also the position of each patch in the topology matters: for instance, we have 
seen that in the network topology described by the chain graph - where all nodes, besides the ones at 
the extremes of the chain, have the same degree - the local dynamics is also influenced by the dynamics 
of the adjacent patches in the graph. Therefore, in hypothetical habitats where there exist many patches 
connected in a linear way, our results suggest that the length of the chain might have a negative role in 
the establishment and in the maintenance of local dynamics. 

Considering the feature of colonization (Section I3.2I ). we have evidenced that, in most network 
topologies, the lack of colonization can be due to the delay of migrating predators with respect to the 
(uncontrolled) local growth of prey, which then leads to the extinction of preys and the prevention of 
the LV dynamics. To effectively measure how strong is the power of the delay, it would be interesting 
to understand whether the local growth of preys can be controlled by inducing their death and thus po- 
tentially allowing the establishment of oscillations. Besides this aspect deserving further investigations, 
our analysis have evidenced that the colonization of empty patches occurs more easily in those patches 
that are adjacent to the patch(es) initialized with the LV complete model. Once more, this highlights the 
relevance of the position of the patch(es) where standard oscillations in preys and predators are already 
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settled at the beginning of the simulation. Indeed, the power of colonization is stronger in the circular 
and complete networks - where the position of the LV complete patch is irrelevant (as the spread of mi- 
grating individuals throughout the network results uniform), and it is weaker in the star network - where 
the position of the LV complete patch is of primary importance (as the spread of migrating individuals 
throughout the network strongly depends on whether the patch is placed at the center or at the tips of the 
star). 

In addition to the investigations that we have presented in this work, further types of analysis that 
we plan to perform on metapopulation systems concern, for instance, the study of the aspects considered 
in this paper (migration, colonization, network topologies, etc.) by assuming other local and global 
dynamics, e.g., the population growth according to the logistic function. Moreover, an interesting issue 
that might be investigated is the synchronization of local population dynamics (e.g. by considering the 
establishment and decay of oscillations in preys and predators) during migration through a given network 
topology, or in the process of colonization. 

Concerning the use of graphs, other relevant questions regard the analysis of the dynamics with 
respect to graph properties, such as different measures of habitat connectivity (centrality indexes) lfl3l 
l26l . In this context, for example, the star graph can resemble the notion of hub (a node with high degree) 
in a typical scale-free network, a structure that is known to be robust to random disturbances but highly 
vulnerable to deliberate attacks on the hubs 021 l3l. 

Another topic of interest concerns the fact that various populations can coexist in a common habitat, 
but have distinct (inter)species dynamics or different dispersal capabilities in that habitat [9]. In cases 
like this, it would be interesting to construct and analyze different metapopulation models, one for each 
target species, according to both the patch-to-patch connections and to the specific population dynamics. 
By comparing and intersecting the results obtained on the distinct network topologies of the common 
habitat derived in this way, it would be possible to determine the locations of the habitat that are most 
important for each species, and thus aid the design of natural reserve systems where we can have the 
most appropriate solution for all species in terms of the maximal improvement of dispersal (reduction 
of species isolation) and the minimal spread of disturbances (diseases, pathogens, invasive species, etc.) 
(361. 

We believe that our modeling approach opens interesting perspectives and can represent an useful tool 
for the investigation of a wide range of properties in metapopulation systems. We expect that applications 
of this model to real cases - characterized by complex habitat networks (where each patch possesses its 
own features of quality, occupancy, connectivity) and different population dynamics - will aid in the 
achievement of important results and new perspective in Ecology. 
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